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Abstract 

Time-dependent stochastic inversion (TDSI) was recently developed for acoustic travel-time 
tomography of the atmosphere. This type of tomography allows reconstruction of temperature 
and wind-velocity fields given the location of sound sources and receivers and the travel times 
between all source-receiver pairs. The quality of reconstruction provided by TDSI depends on 
the geometry of the transducer array. However, TDSI has not been studied for the geometry 
with reciprocal sound transmission. This paper is focused on three aspects of TDSI. First, the 
use of TDSI in reciprocal sound transmission arrays is studied in numerical and physical 
experiments. Second, efficiency of time-dependent and ordinary stochastic inversion (SI) 
algorithms is studied in numerical experiments. Third, a new model of noise in the input data 
for TDSI is developed that accounts for systematic errors in transducer positions. It is shown 
that (i) a separation of the travel times into temperature and wind-velocity components in 
tomography with reciprocal transmission does not improve the reconstruction, (ii) TDSI yields 
a better reconstruction than SI and (iii) the developed model of noise yields an accurate 
reconstruction of turbulent fields and estimation of errors in the reconstruction. 

Keywords: travel-time tomography, inverse problems, acoustic imaging, remote sensing, 
tomography of the atmosphere, stochastic inversion 


1. Introduction 

Acoustic travel-time tomography of the atmosphere allows one 
to reconstruct temperature and wind-velocity fields in a certain 
plane or volume given the coordinates of sound sources and 
receivers and travel times of sound propagation between all 
source-receiver pairs [1-6]. In outdoor acoustic tomography 
experiments [7], the number of sources and receivers has been 
relatively small (less than 20) so that the number of measured 
travel times has also been small. Therefore, one of the main 
issues in acoustic tomography of the atmosphere, similarly to 

5 Present address: Signature Physics Branch, US Army ERDC-CRREL, 72 
Lyme Rd, Hanover, NH 03755-1290, USA. 


analogous problems in other areas [8-10], is the formulation 
of robust and accurate inverse algorithms for reconstruction 
of the temperature and wind-velocity fields from a limited 
amount of data (the travel times and transducers’ coordinates). 

To address this issue, a time-dependent stochastic 
inversion (TDSI) algorithm was recently developed for travel¬ 
time tomography of the atmosphere [11-14] . The main idea of 
TDSI is to measure the travel times repeatedly (i.e., at different 
times) and to use all these measurements for reconstruction 
of temperature and wind-velocity fluctuations at an arbitrary 
moment of time. This can be accomplished by taking into 
account temporal covariances of the fluctuations. TDSI takes 
into account both spatial and temporal covariances and, hence, 
is an improvement over the ordinary stochastic inverse (SI) 
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that accounts for spatial covariances only [1, 15]. For the 
problem considered in this paper, significant correlations exist 
in the atmospheric surface layer over time scales of tens of 
seconds and spatial scales of tens to hundreds of meters. TDSI 
might also be considered for tomographic applications on a 
larger scale, such as global atmospheric tomography using 
infrasound. In this case, it could take advantage of correlations 
involving synoptic weather patterns, which have correlation 
time scales of days and spatial scales of hundreds of kilometers. 

The general idea of TDSI is successfully applied in other 
areas. For example, a similar approach is used in satellite 
altimetry of the ocean surface to interpolate measurements 
along satellite tracks to other spatial and temporal points 
[16-18]. In image processing, this approach is used to improve 
the image quality from several blurred or noisy frames [ 19, 20] . 
A similar algorithm is used in medical tomography where 
it is referred to as the vector Wiener filter [21]. TDSI is 
also somewhat similar to the Kalman filter [22, 23] which is 
applied in many areas including tomography [24-27]. The 
distinction is that TDSI is not recursive (i.e., yields a closed- 
form reconstruction), does not make assumptions regarding 
stochastic variables that adopted in Kalman’s theory (e.g., 
does not require formulation of a problem in state-space terms 
which could be a difficult task) and uses data obtained at 
arbitrary times for the reconstruction of fluctuations at other 
arbitrary time (in this paper, for example, the previous, current 
and future measurements relative to the time of reconstruction 
are used). 

The question about what functions can be used to describe 
spatial covariances of isotropic, homogeneous and statistically 
stationary atmospheric turbulence, which is considered in this 
paper, is well studied in the literature (e.g., [28-32]). The 
spatial-temporal covariance functions can be obtained from 
spatial ones with the use of frozen or locally frozen turbulence 
hypotheses [11, 12, 28]. However, these functions may not 
coincide with spatial-temporal covariance functions of actual 
turbulence since (i) the theoretical and actual spatial functions 
may differ, and (ii) the frozen or locally frozen hypotheses 
may not describe actual turbulence adequately. As a result, the 
TDSI algorithm may not yield more accurate reconstruction 
than an ordinary stochastic inversion algorithm, which utilizes 
only spatial covariance functions. The advantage of TDSI in 
comparison with SI was demonstrated in [11] for numerically 
created rigidly frozen turbulence. However, it has not been 
shown that TDSI provides a more accurate reconstruction in 
the case of realistic turbulence which may or may not be frozen. 

The goals of the paper are threefold. First, 
application of TDSI to acoustic tomography with reciprocal 
sound transmissions is studied. Reciprocal transmission 
tomographic arrays have been used successfully in underwater 
and atmospheric tomography for many years [33-36]. There 
are two advantages of such arrays. First, they allow one 
to increase the number of measured travel times for a 
given number of transducer stations since the total number 
of sound sources and receivers increases. Second, they 
allow one to separate the originally measured travel times 
into two components: travel times due to temperature and 
travel times due to velocity. This makes the contribution 


of weaker fluctuations to the data apparent. Therefore, one 
can solve two problems (one for the reconstruction of the 
temperature field and another one for the reconstruction of 
the wind-velocity field) more effectively than the original 
problem, in which temperature and wind-velocity fields are 
reconstructed simultaneously. However, such a separation 
reduces the number of travel times for each of the fields 
by a factor of two. Therefore, the analytical advantage 
of independent field reconstruction may not lead to better 
quality. For some algorithms, such an analytical separation 
led to a better reconstruction (or was a presupposition of 
some specific algorithms) [6, 36, 37]. However, it has 
not been studied whether such an approach improves the 
reconstruction with the use of stochastic algorithms. In this 
paper, reciprocal transmission arrays are studied for time- 
dependent and ordinary stochastic inversion algorithms. 

The second goal of this paper is to verify that 
TDSI yields a better reconstruction than SI for realistic 
turbulence. For this purpose, time-dependent, horizontal 
slices through temperature and wind-velocity fields from a 
high-resolution large-eddy simulation (LES) of an unstably 
stratified atmospheric boundary layer were used as a surrogate 
atmosphere. The LES solves an approximate, filtered form 
of the nonlinear fluid dynamic (Navier-Stokes) equations that 
explicitly describes the large-scale structure of the flow. More 
detailed information about the LES used in this paper can be 
found in [38]. 

Finally, in this paper, a more complete model of noise in 
the input data for SI and TDSI is developed. This model 
assumes that the errors (noise) in the transducer positions 
is systematic; i.e. it remains the same for measurements at 
different times for a given ray path. Also, the errors in the 
mean field reconstruction affect data corresponding to all travel 
paths. These improvements result in a non-diagonal noise 
matrix, in contrast to earlier publications [11-14]. 

All derivations are implemented for a two-dimensional 
problem where all sources and receivers lie in one horizontal 
plane, and two-dimensional temperature and wind-velocity 
fields are subject to reconstruction. The derived results remain 
valid for the three-dimensional case. 

The paper is organized as follows. The separation of 
travel times into temperature and wind-velocity travel times 
is presented in section 2. The basic equations for the 
reconstruction of temperature and wind-velocity fields are 
given in section 3. The numerical experiment with the use 
of LES fields is described in section 4. Then, application of 
TDSI to real experimental data from a reciprocal transmission 
array is demonstrated in section 5. Summary and conclusions 
are presented in section 6. 

2. Separation of travel times 

In this section, it is shown how the temperature and wind- 
velocity contributions to the travel times can be separated. 
The linearized equation for travel times in two-dimensional 
acoustic tomography of the atmosphere in the horizontal (x , y ) 
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plane is given by the following equation [11]: 

. _ U L _ u 0 (t)s ix + v 0 (t)s iy \ 

'■ u _ co(0 v Co(0 / 

- At f dl I Arf T (r, t) + u(r, t)s ix + v(r, t)s iy |, 
CqU) J/,, 1 27 0 (f) J 

(1) 

where tf is the travel time of a sound propagating along path 
Li, i = 1,2,...,/,/ being the total number of travel paths 
(if S is the number of sound sources and R is the number of 
receivers, then I = SR), Vo = uq£ x + VQe y is the vector of 
wind velocity spatially averaged over a tomographic plane, e x 
and are the basis vectors of a Cartesian coordinate system, 
Co and Tq are adiabatic sound speed and temperature averaged 
over a tomographic plane, V = ue x + ve y and T are fields 
of fluctuations, s/ = Si x e x + Si y e y is a unit vector in the 
direction of the ith ray propagation, r = xe x + ye y is a spatial 
vector and t is time. The integration is taken along path L t 
which, in this approximation, is a straight line. The full fields 
of temperature and wind velocity in the tomographic plane, 
T (r, t) and V(r, t), can be found as the sum of their spatial 
mean values co(0 and Vo (t) and fluctuation fields T(r,t) and 
V(r, t): 

T (r, t) = 7b(0 + T (r, t ), V(r, t) = V 0 (f) + V(r, t). (2) 

The mean values Co and Tq are connected by the following 
relationship: 

Cq(0 = yR a T 0 (t), (3) 

where y = c p /cy ~ 1.41 is the ratio of specific heats c^ 
and Cy , and is the gas constant for dry air. Note that T 
is the acoustic virtual temperature which relates to ordinary 
thermodynamic temperature T th as T = T th (1 + 0.511C), 
where C is the concentration of water in air [31]. Note that, 
for small C, C = q, where q is the specific air humidity. The 
temporal variations of the mean fields could be significantly 
greater than the magnitude of the fluctuations. For example, in 
the field experiment described in [12], Tq varied from 15.9 to 
16.4 °C, uq from —0.3 to —1.3 m s -1 , and vq from —1.6 to 
—2.1 m s -1 during 10 min while the corresponding standard 
deviations of the fluctuations were a T = 0.27 °C and 
oy — 0.28 m s -1 . 

The goal of acoustic travel-time tomography is to 
reconstruct fields T 0 (t), V 0 (f), T(r, t) and V(r, t) given the 
coordinates of sound sources and receivers (which determine 
Sf) and travel times tf(t) for all rays. 

Note that there may exist such wind-velocity fields V 
that are ‘invisible’ for travel-time tomography in the sense 
that they do not contribute to the travel times tf [34, 35]. 
To detect these fields some additional measurements, besides 
travel times, are suggested [34]. In this paper, no additional 
measurements are used for reconstruction of wind velocity 
with the use of stochastic inversion algorithms. The results of 
such a reconstruction are discussed in section 5. 

A reciprocal ray setup provides the travel times for each 
pair of rays simultaneously propagating in opposite directions. 
This allows an analytical separation of the measured travel 
times on two subsets: the travel times tf T , which are affected 


only by T, and travel times tf v , which are affected only 
by V. Let n and k be indices for reciprocal rays: n = 
1, 2,..., N,k = 1, 2,..., N, and N = 1/2. Their lengths 
are the same, L n = L^, but their directions are opposite: 
s k = —s n . Taking the sum and difference of the measured 
travel times tf and tf and using equation (1), one obtains 


CAO = 


+ L n 


1 


Co(0 2c 0 (t)T 0 (t) 


L 


T(r,t)dl, 

( 4 ) 


*tr /*\ £(*)-£(*) 

l nV\ l > ~ 2 

_ L n (u 0 (t)s 

nx "I" Vo(f)Sny) 

Cq(0 

+ f dl(u(r,t)s nx + v(r,t)s ny ). (5) 

clit) J Ln 

The advantage of equations (4) and (5) in comparison with 
the original equation (1) is that they allow independent 
reconstruction of temperature and wind-velocity fields. 
However, the amount of data for such reconstruction is 
two times less. Therefore, it is worthwhile comparing the 
reconstructions based on equation (1) and on equations (4) 
and (5). 


3. Reconstruction of the mean fields 

The reconstruction of mean fields uo(t),vo(t) and To(t) 
will be based on the same formalism that was utilized in 
[11-13]. However, the direct application of the formalism, 
developed in these references, is not possible since the starting 
equations (i.e., equations (4) and (5)) differ from the original 
one (equation (1)). 

The first step is to neglect fluctuations in equations (4) and 
(5). This results in omitting the integrals in these equations. 
Physically, this means that one reconstructs the uniform fields 
co(0 and Vo(0 which match in the best way the known travel 
times tf T (t) and tf v (t). The value of Tq is calculated from cq 
with the use of equation (3). The second step is to use travel 
times at time t to reconstruct the mean fields at the same time 
t. Then, time ns a parameter and will be omitted until the end 
of this section. The fact that travel times depend on t will be 
used in the following section in which fluctuations T(r,t) and 
V(r, t) will be reconstructed by TDSI. 

3.1. Reconstruction of mean temperature 

After omitting the integral in equation (4), the equation for cq 
can be rewritten in the following form: 

CQ = d nT , (6) 

where d n j = L n /tf T . The well-known least-squares solution 
Co and the estimation of its standard deviation a Co for such a 
problem are given by the following formulae [15]: 

(7 ) 
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Vm — 


'Yln= 1 WnT ~ Q)) 2 

N(N ~ 1) 


( 8 ) 


The reconstruction of mean temperature To is obtained 
from co with the use of equation (3). 


The errors <xj? and cf? o take into account the errors in travel 
time measurements, transducer positions, and omitting the 
integral term in equations (4) and (5). (Note that errors in 
travel-time measurements refer to a specific technique of the 
measurement; they are independent of the errors caused by the 
transducer position uncertainty or omitting the integral.) 


3.2. Reconstruction of mean wind velocity 


Once the value of Co is reconstructed, one can use equation (5) 
to reconstruct the mean values of wind-velocity components. 
Neglecting the fluctuations in equation (5), one can formulate 
the problem in matrix notation as follows: 


where 


and 


Gf = dy, 

$lx ^ly 

G= : : , 

_$Nx SNy_ 

f = [m 0 ; v 0 ], 

dy = ^o\flv/ ^1’ * • • ’ i'NV / 


( 9 ) 


GO) 

( 11 ) 

( 12 ) 


The semicolon between elements indicates that these elements 
are arranged in a column so that f and dy are column vectors. 
In equation (9), the matrix G and the vector d v are known, 
while the vector f is subject to estimation. Usually, A > 2, and 
the problem seems to be overdetermined since the number of 
equations is greater than the number of unknowns. However, 
the matrix G could be ill-conditioned due to an unfortunate 
location of the transmitters and receivers. In this case, the 
straightforward least-squares solution of equation (9), which 
implies inversion of matrix (G r G), cannot be employed 
accurately. Therefore, it is worthwhile using the pseudo¬ 
inverse [39] matrix G” 1 , in terms of which the solution of 
equation (9) is given by the following formula: 


f = G; 1 dy. 


( 13 ) 


If the G matrix is well conditioned (i.e., its rank equals 2), 
then G ” 1 coincides with the least-squares estimator: G ” 1 = 
(G r G) - 1 G r , and equation (13) yields the least-squares 
solution of equation (9). If the G matrix is ill conditioned 
(i.e., its inverse condition number equals or approaches 
zero), then it corresponds to the situation when one of the 
unknowns («o or vq) is overdetermined while the other is 
underdetermined. In this case, G ” 1 yields the least-squares 
solution for the overdetermined unknown and the minimal 
second norm solution for the underdetermined one. 

Similar to equation ( 8 ), the mean squared errors a] = 
[o^ o , of the solution f can be estimated with the help of 
the residual estimation [ 12 ]: 

??=4diag(G; 1 [G; 1 ] 7 '). (14) 

where s d is the estimated variance in the data dy: 

2 _ (dy - Gf) r (dy - Gf) 


4. Reconstruction of the fluctuations 

In this section, the fluctuations T(r,t) and V(r, t) will be 
reconstructed with the use of the TDSI algorithm. The 
equations of this algorithm for travel-time tomography of 
the atmosphere were derived in [11]. Here, it is shown 
how to apply TDSI to reciprocal transmission transducer 
arrays, when equation ( 1 ) can be split into two equations, 
equations (4) and (5). Application of TDSI to acoustic 
tomography with reciprocal sound transmission is similar 
to that for non-reciprocal sound propagation developed in 
[11-14]. Therefore, it is worthwhile presenting a brief 
overview of the main results obtained in these papers. 


4.1. Time-dependent stochastic inversion 


The starting equation for the reconstruction of the fluctuations 
T( r, t) and V(r, t) is obtained from equation (1) where the 
first term is calculated with the use of the reconstructed values 
< co(t),uo(t ) and To (0- The remaining part can be written in the 
following form: 

q,(t) = qoi(t) + ^(t), (16) 


where 


q;(t) = Li(c 0 (t ) - s, • %«)) (17) 


got (t) is noise-free data: 


<7o i(t) = 


L 


~ cp(t) 

.27b(f) 


T (r, t) - Si • V(r, t) 


dZ, 


( 18 ) 


and & (t) represents noise in q t (t) due to the errors in travel 
time measurements, transducer position measurements and the 
reconstruction of mean fields. 

Suppose that the travel times t tr (0 = [tf(t); ...; tf(t)] 
were measured at times t\, h ,... ,tg. For each 4 (k = 
1 , 2 , ..., 2 ), the mean fields are estimated similarly to the 
technique described in section 3 (see also [11, 12, 14]), and 
vectors q( 4 ) are formed. The key idea of TDSI is to use all Q 
vectors q(4) to reconstruct fluctuations T( r, 4 ) and V(r, 4 ) 
at arbitrary time 4 . 

Let m(4) be a vector of models which are subject to 
reconstruction: 


m(4) = \T( 14, 4 ); ...; T(rj, 4 ); u(r u 4 ); ...; u(r j9 4 ); 

v(r u to); v(rj,t 0 )], (19) 

where J is the total number of spatial points within 
a tomographic plane where the fluctuations should be 
reconstructed and d is a vector of data for this reconstruction: 


d = [q(4);qfe);-..;q(^)]. ( 20 ) 
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Then, the optimal stochastic estimation m (to) of models m(to) 
is given by the following equation [ 11 ]: 

m(fo) = R 111 C |R dd 'd, (21) 

where R md = (md 7 ) and R dd = (dd 7 ) are model-data 
and data covariance matrices (the angular brackets () denote 
the mathematical expectation). The estimation m(^ 0 ) given 
by equation ( 21 ) minimizes the expected squared errors (e 2 ) 
(j = 1 , ...,/) for each spatial point of the reconstruction 
[11,15 \.(t))=((m j -m j ) 2 ). 

The elements of R m d and Rdd matrices are calculated 
with the use of spatial-temporal covariance functions of 
temperature and wind-velocity fluctuations B TT (r', t'; r", t") 
and Bij(r\ t'\ r", t") (here, i and j indices indicate two 
components of wind-velocity fluctuations V(r, t)). Formulae 
for the calculation of the elements of noise-free matrices 
Rmdo an d Rd 0 d 0 are given by equations (18)—(21) from [ 11 ]. 
For convenience, these formulae are presented in appendix 
A. Matrices R m d and Rdd are obtained from R m d 0 and 
Rd 0 d 0 by taking into account noise in the input data, as 
shown in appendix A. To find analytical formulae for the 
spatial-temporal covariance functions B TT (r', t'\ r ", t") and 
Bij(r', t'; r", t"), the hypothesis of frozen or locally frozen 
turbulence can be used [11-13, 28, 29]. These formulae can 
be found in appendix B. 

The mean squared errors of m (to), which are the errors of 
the reconstruction of m(f 0 ), are given by the diagonal elements 
of the error covariance matrix R e€ [ 11 ]: 

Ree = Rmm — Rinded ^md’ (22) 

where R mm = (m(to)m T (to)) is the model covariance matrix. 

4.2. Error analysis in the data for TDSI 

In this subsection, the covariance matrix of noise in the input 
data for TDSI, R^, is calculated. Knowledge about noise 
in the input data is essential for TDSI since it appears in the 
Rdd matrix, as shown in appendix A. The elements of the 
R^ matrix are given by (§i(O£r(* 0 )- To find these elements, 
one can take variations of equation (17) for different rays i 
and i' at different times t and t\ multiply them, and take 
the mathematical expectation. Taking into account that for 
straight rays LiS ix = (x Ri - x Ti ) and LiS iy = (y Ri - yn), 
where (. x Ri ,y Ri ) and (x Ti , y Ti ) are the coordinates of a receiver 
and a transmitter of the ith ray, correspondingly, one has: 

£/(0 = &qi(t) = Li8co(t) + cb(05Lj -u 0 (t)(8x Ri - 8x n ) 
-Do (t)(Sy Ri - 8y n ) - (x Ri - x Ti )8u 0 (t) 

- (yRi ~ yn)Svo(t) ~^l(t)8tf(t) - 2co(t)t*(t)8c 0 (t). 

(23) 

A similar equation can be written for a different path i' at 
a different time t'. Since L 2 = (x Ri — xn) 2 + (y R i — yn) 2 , 
then: 

« T ( x Ri ~ x n) (8x Ri — 8x T i) + (yRi — yn) (8yRi — 8yn) 

8Et = --- 

E n 

= s ix (8x Ri - 8x Ti ) + s iy (8y Ri - 8y n ). (24) 


Taking into account that Do tf ~ L t , many of the cross¬ 
correlations equal zero, and neglecting the terms of order 
(Vo/Do) 2 , one has 

(^i(t)^(t')) & 2 a 2 (D 0 (ODo(D) - lc 0 (t)uo(t') +'co(t')uo(t)]s ix 
- [D 0 (ODo(D) + c‘o(t f )vo(t)]siy)8 w 

+ LiLy (&co ^ix^i'x^uo SiySi'y&vo)^tt' + &ii'&tt'i 

(25) 

where a 2 is the variance of the errors in transducer positions 
(a 2 = ((8x Ri ) 2 ) = {(8x n ) 2 ) = (( 8y Ri ) 2 ) = ((5y r/ ) 2 )), 
'a 2 cr 2 and cr 2 are variances of errors in the reconstruction 

Co ’ u o’ Vo 

of Co, uo, and Vo, correspondingly (see section 3 and similar 
formulae in [ 11 , 12 ]), cr 2 is the variance of errors in the travel 
time measurements, and 8 W and 8 tt > are Kronecker’s delta 
symbols. While deriving equation (25), it is assumed that 
uncertainties in the transducer positions are systematic in the 
sense that they do not change with time for any given path. 
As a result, the noise for the same paths is correlated for data 
measured at different times t and t' (the first term in equation 
(25)). Also, the errors in the reconstruction of mean fields 
at a given time t affect all data measured at this time (the 
second term in equation (25)). These two features result in 
the non-diagonal noise matrix R^ which should be added to 
the noise-free data covariance matrix Rd 0 d 0 ( see appendix A). 
Finally, in the particular case t = t' and i = i', equation (25) 
yields more accurate estimation of the variance than that 
presented in [ 12 ]. 

4.3. TDSI for reciprocal arrays 

The starting equations for the reconstruction of temperature 
and wind-velocity fluctuations are derived from equations (4) 
and (5): 

qnT (t) = qOnT (t) + %nT (t), 4nV (0 = <10n V (0 + %nV (0> (26) 

where q nT and q nV are calculated with the use of the 
reconstructed mean values as described in section 3: 


q nT (t) = L„co(t) -Co(t)t* T (t), (27) 

q n v(t) = Co(OCv(0 “ L n(UoSnx +^ny), (28) 


the noise-free data qo nT and qo nV are given by the following 
expressions: 


qonr(t) = 


-L 


/JI 


Do ( 0 _ 

(0 


T( r, t) 


d /, 


qonv(t) = / (u( r, t)s nx + v(r, t)s ny ) d/, 


(29) 

(30) 


and % nT and § nV represent the noise in the data due to the errors 
in measurements of travel times, positions of the transducers 
and reconstructions of the mean fields. The variances of these 
noises can be calculated similarly to the derivation shown in 
the previous subsection, but the analysis should be based on 
equations (27) and (28) instead of equation (17): 

(^nT(t)^n'T(t')) ~ 2& 2 Co(t)Co(t')8 nn ' 

+ L n L n >af o S tt : + \v%(t)a?8 nn '8„', ( 31 ) 
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(b) 



0 Sources 
□ Receivers 
Rays 


Figure 1 . Tomographic arrays: ( a ) for the numerical experiment and ( b ) for the physical outdoor experiment. 


(£nv(0§n'v(f')) ~ 2a^(Uo(tyUo(t') + %(t)Vo(t'))8 n n' + L n L n ’ 

x {^nx^n'x^uQ + SnySn'yGy 0 ) &tt' + 2^0(0^ &nn'&tt'i (32) 

where indices n and ft' denote two reciprocal paths. Note that 
the variance of (the last term in equation (32)) is one-half 
of the variance of measured travel times. This is also true for 

l nV' 

Comparing equation (18) with equations (29) and (30), 
one concludes that the general TDSI algorithm can be used 
for reconstruction of either the temperature or the wind- 
velocity fields with few modifications. For the temperature 
field T, one should use a different data vector d T = 
[qr( 0 ); qrfe); ...; q t^q)], where q^ are given by equation 
(27), set the wind-velocity fluctuations to zero (cry = 0 in 
formulae (A.3) and (A.4) for the calculation of R m( j 0 and 
Rdodo matrices in appendix A) and use equation (31) instead of 
equation (25) to calculate the error matrix, which is needed in 
order to find the matrix Rdd- Correspondingly, to reconstruct 
the field of wind-velocity fluctuationsV, one should use the 
data vector d v = [qy (L); qvfe); ...; qv(^)], where q v are 
given by equation (28), set cr T = 0 in formulae (A.3) and 
(A.4), and use equation (32) to calculate the noise. 

5. Numerical experiment 

In this section, a numerical two-dimensional tomography 
experiment is described. The primary goal of this numerical 
experiment is to study whether the analytical separation 
of travel times in reciprocal transmission arrays, given 
by equations (4) and (5), improves the reconstruction in 
comparison with the general TDSI algorithm based on 
equation (1). The secondary goal is to compare the 
reconstruction quality of TDSI and ordinary SI algorithms 
on non-frozen turbulence. 

For these purposes, the original temperature and wind- 
velocity LES fields were used. The particular simulation 
analyzed here involves unstable atmospheric stratification. 


Five time frames ( Q = 5) of two-dimensional temperature 
and wind-velocity fields were employed with a time interval of 
4.6 s at height z = 13.75 m. The spatial resolution of LES was 
4 m in the horizontal plane. More detailed information about 
LES fields used in the numerical experiment described can be 
found in [38]. Using these LES fields, the travel times were 
calculated for each time frame in accordance with equation ( 1 ) 
for the reciprocal transmission array shown in figure 1(a). The 
array consisted of eight sound sources (SI,..., S 8 ) and eight 
receivers (R1,..., R 8 ) arranged along the perimeter of a square 
with an 80 m length side. (Note that a similar array is currently 
under construction at the National Oceanic and Atmospheric 
Administration, Boulder, CO.) The array provided 56 travel 
times (eight sources x eight receivers — eight paths of zero 
length) for each time frame (28 reciprocal travel paths). To 
take into account the noise in travel-time measurements, the 
normally distributed white noise with a t =5 /zs was added 
to the calculated travel times (for each path at each time 
frame). This magnitude of noise corresponds to 3% of 
errors in the data for the reconstruction of fluctuations (see 
equation (16)). Moreover, during the calculations of travel 
times, the transducer positions were randomly distorted by 
normally distributed white noise with a r = 0.01 m to represent 
the uncertainty in the transducer positions. These random 
deviations from the original locations remained the same for 
all five time frames to imitate the systematic errors. The 
objective of the experiment was to reconstruct the temperature 
and wind-velocity fields at the third time frame (to = 3). 

The reconstruction of mean fields has been implemented 
for the separated travel times as described in section 3, and the 
reconstruction based on the travel times without separation 
is described in [11-13]. The true and reconstructed mean 
fields are presented in table 1 , where ‘modified alg.’ and 
‘general alg.’ refer to the reconstruction with and without 
travel-time separation, correspondingly. According to table 1 , 
the values of the mean fields reconstructed by both modified 
and general algorithms are the same. The only difference is 
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Figure 2. Original LES and reconstructed temperature and wind-velocity fields at frame 3 in the numerical experiment, (a), (d) and (g): 
original fields of temperature and two components of wind velocity. ( b ), ( e ) and (h): fields reconstructed by TDSI. (c), (/) and (/): fields 
reconstructed by SI. 


Table 1 . Actual and reconstructed mean temperature and 
wind-velocity fields. 


Fields 

To (°C) 

u 0 (m s l ) 

v 0 (ms l ) 

True 

28.55 

0.80 

1.82 

General alg. 

28.56 ± 0.07 

0.80 ± 0.05 

1.82 ± 0.05 

Modified alg. 

28.56 ± 0.01 

0.80 ± 0.08 

1.82 ± 0.08 


in the estimated errors of the reconstruction. But these errors 
are still of the same order. Another conclusion from table 1 
is that the reconstruction of the mean fields is very accurate. 
Comparing the true and reconstructed values, one can see that 
the actual discrepancy is of order 0.01 °C for temperature and 
less than 0.01 ms -1 for wind velocity. 

The reconstruction of the fluctuations was implemented 
by TDSI and SI algorithms with and without analytical 
separation of travel times. For the TDSI, the data from all five 
time frames were used to reconstruct the fields at time to, while 
SI utilizes the data measured at the same time frame ^o- The 


variances and correlation lengths of turbulence for stochastic 
inversion algorithms were estimated from the original LES 
fields: oj = 0.08 °C, <r v = 0-38 m s -1 and l T — l v — 14 m. 
Figure 2 shows the original and full fields T , u, and v 
reconstructed without travel-time separation. The original 
LES fields are presented in figures 2(a), (d) and (g). The fields 
reconstructed by the TDSI algorithm are shown in figures 2(b), 
(e) and (h). The reconstruction with the use of SI is presented 
in figures 2(c), (/) and (i). Comparing these figures, one 
notes a remarkable improvement in the reconstruction done 
by TDSI in comparison with that done by SI. In the case 
of SI, the temperature field is reconstructed rather poorly. 
Neither the shape of turbulent eddies nor their magnitude is 
reconstructed sufficiently well (see figure 2(c)). In contrast, 
the TDSI algorithm (figure 2(b)) captures both the shape and 
the magnitude even though the smallest details are still omitted. 
The SI algorithm reconstructs the wind-velocity components 
better than the temperature field (see figures 2(f) and (/)), but 
still not as well as the TDSI algorithm (figures 2(c) and (h)). 
For example, the shape of blue eddies in the u and v fields is 
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Table 2. Actual and expected RMSE of the reconstruction by TDSI 
and SI algorithms. 


RMSE 

T (°C) 

u (ms l ) 

v (ms x ) 

TDSI General algorithm 

Actual 

0.05 

0.233 

0.215 

Expected 

0.08 

0.235 

0.208 

Modified algorithm 

Actual 

0.05 

0.226 

0.209 

Expected 

0.05 

0.233 

0.207 

SI General algorithm 

Actual 

0.07 

0.313 

0.286 

Expected 

0.10 

0.330 

0.267 

Modified algorithm 

Actual 

0.08 

0.307 

0.284 

Expected 

0.08 

0.326 

0.266 


not reconstructed correctly by SI while TDSI reconstructs it 
rather well. 

These observations are summarized in table 2, which 
represents actual and expected root mean squared errors 
(RMSE) for the general (without the separation of travel times) 
and modified (with travel times being separated) algorithms 
for the case of SI and TDSI. The actual RMSE is calculated 
with the use of an actual discrepancy between the original and 
reconstructed fields. The expected spatially averaged RMSE 
of the full fields are calculated with the use of equations (8) and 
(14) for the mean fields and equation (22) for the fluctuations. 
Comparing the actual errors of the general and modified 
algorithms (lines 1 and 3 for TDSI and lines 5 and 7 for SI 
in table 2), one concludes that these algorithms yield almost 
identical results. This means that the analytical separation of 
the travel times does not improve the reconstruction quality 
of the stochastic inversion algorithms. On the other hand, if 
one compares the actual errors of TDSI and SI reconstructions 
(lines 1 and 5 for the general algorithm and lines 3 and 7 
for the modified one), then it is noticeable that TDSI always 
outperforms SI. This is clearly seen in figure 2 as well. Finally, 
the comparison of the actual and expected errors for the 
same algorithms reveals that the technique used for the error 
estimation yields accurate and reliable results. 

It is interesting to point out that the ‘invisible’ wind- 
velocity fields, which are mentioned in section 2 and [34, 35], 
did not introduce significant errors in the reconstruction. The 
actual RMSE are small and in a good agreement with the 
expected RMSE. Two explanations are plausible. First, these 
‘invisible’ fields were not reconstructed but they were much 
weaker than ‘visible’ ones, which were reconstructed from 
travel times, so that they could not distort the reconstructed 
fields significantly (within the actual RMSE of reconstruction). 
This point can be supported by the following reasoning. The 
use of Taylor’s or locally frozen turbulence hypotheses in 
TDSI may be viewed as effective moving the transducer 
array through motionless turbulence. This is equivalent to 
having some transducers not only along the borders of the 
tomographic area but also inside. The data obtained from 
such an extended transducer array at several time frames 
may reveal ‘invisible’ fields that exist longer than the time 
interval between the frames. Since the time interval between 


Table 3. Mean temperature and wind-velocity fields reconstructed 
in the outdoor experiment. 


Fields 

To (°C) 

w 0 (m s l ) 

v 0 (m s l ) 

General alg. 

27.30 ± 0.06 

-0.15 ±0.05 

0.27 ± 0.05 

Modified alg. 

27.30 ± 0.07 

-0.15 ±0.05 

0.27 ± 0.05 


Table 4. Expected RMSE of the reconstruction by TDSI algorithm 
in the outdoor experiment. 


RMSE 

T (°C) 

u (ms l ) 

v (ms l ) 

General alg. 

0.08 

0.26 

0.25 

Modified alg. 

0.09 

0.26 

0.25 


two consecutive measurements can be quite small (it is an 
adjustable parameter in experiments), only really short-time 
and weak ‘invisible’ fields may remain unrevealed. Second, 
the ‘invisible’ fields may be reconstructed by stochastic 
inversion algorithms (in full or partially). As shown in 
section 4, these algorithms incorporate additional statistical 
information about true fields. This information supplements 
the travel-time measurements (in analogy to the additional 
measurements proposed in [34]). Then, the ‘invisible’ fields 
are reconstructed by stochastic inversion algorithms as a 
complement to the ‘visible’ fields so that the total fields would 
have specific correlation properties. 

6. Physical outdoor experiment 

This section describes a physical outdoor experiment with 
a reciprocal transmission array carried out as a part 
of the experiment STructure of turbulent transport under 
INHOmogeneous surface conditions (STINHO) on June 17, 
2002, in Lindenberg, Germany [7]. The array consisted of 
eight sound sources and eight sound receivers arranged along 
the perimeter of 250 m x 300 m as shown in figure 1(b). The 
travel times were measured each minute. The order of errors 
in the measurements was o t = 0.1 ms for the travel times 
and a r — 0.01 m for the transducer positions. The fields at 
1753 UTC (universal time coordinated) were the subject of 
reconstruction. 

The reconstruction was implemented in the same manner 
as it was for the numerical experiment except that there 
were no true fields and actual errors of the reconstruction. 
However, the results of the numerical experiment suggest that 
the reconstruction by TDSI is quite accurate, and the expected 
errors yield an accurate estimation of the actual ones. 

The reconstruction of mean fields by general and modified 
algorithms is presented in table 3. Similar to the numerical 
experiment, the reconstructed values and expected RMSE are 
the same for the general and modified algorithms. As one 
can see, the expected RMSE are 0.06 °C for temperature and 
0.05 m s -1 for wind velocity. This accuracy is sufficient for 
most meteorological applications. 

For the TDSI algorithm, the variance of temperature 
fluctuations was estimated with the use of an in situ sensor 
located within a tomographic plane: a T = 0.06 °C. The 
variance of the wind-velocity fluctuations was not measured 
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Figure 3. Physical outdoor experiment, (a)-(c): temperature and two components of wind-velocity fields reconstructed by TDSI. The black 
lines in (a) indicate the reconstructed wind, (d)-(f): the expected errors of the reconstruction shown in (a)-(c), correspondingly. 


in this experiment but was measured in similar outdoor 
experiments [7, 12]. In the experiment described, the typical 
value ay = 0.3 m s -1 had been used. The correlation lengths 
of the fluctuations, which define a desirable characteristic size 
of the reconstructed turbulent eddies, were l T = l v = 40 m. 
The travel times from three consecutive frames were used 
(1752, 1753 and 1754 UTC) to reconstruct the fluctuations 
at 1753 UTC. The reconstructed full fields (the mean fields 
plus the fluctuations) are presented in figures 3 (a)-(c). The 
black lines in figure 3(a) represent the reconstructed wind. 
The expected RMSE of the reconstruction are shown in 
figures 3 (d)-(f). Note that the errors are smaller than the 
range of the fluctuations, which means that the temperature 
and wind-velocity eddies which are seen in figures 3 (a)-(c) 
are reconstructed reliably. 

7. Summary and conclusions 

The main goal of this paper was to study whether reciprocal 
transmission arrays improve tomographic reconstruction of 
atmospheric temperature and wind-velocity fields when 
implemented by stochastic inversion algorithms. It is known 
that, for the linearized problem, reciprocal transmission 
arrays allow one to separate the measured travel times into 
two components, one of which depends on the temperature 
field only while another depends on the wind-velocity field. 
Also, it is well known that such an approach improves the 
reconstructions by some other algorithms since it increases the 


data/unknowns ratio. However, application of this technique 
to stochastic inversion algorithms has not been studied. The 
results obtained in this paper suggest that the modified and 
general stochastic inversion algorithms (with and without 
analytical separation of travel times, correspondingly) yield 
practically identical reconstructions. This conclusion has been 
verified for the TDSI and SI algorithms on the numerical and 
real experimental data. 

Another new contribution of this paper is an improved 
model of noise in the data for reconstruction of the fluctuations 
with the use of SI and TDSI. The present formulation, unlike 
previous ones, accounts for systematic noise in the transducer 
positions, which better corresponds to a real experiment. As 
a result, the covariance matrix of noise is no longer diagonal. 
This extension provides more accurate reconstruction (and 
estimation of the errors in reconstruction) than previously. 

Finally, the TDSI algorithm has been tested for the first 
time on the original LES fields which were not frozen. The 
TDSI algorithm accounts for the correlation of the temperature 
and wind-velocity fields not only in space but also in time. 
To find the spatial-temporal covariance functions of the 
fluctuations, a hypothesis of locally frozen turbulence was 
utilized. Since this hypothesis makes certain assumptions 
about turbulence which may not occur in real turbulence, the 
TDSI algorithm might not result in a better reconstruction 
than the ordinary SI algorithm, which is free from these 
assumptions. The numerical results obtained in this paper 
(figure 2 and table 2) show that the TDSI reconstruction, which 
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utilizes the locally frozen turbulence hypothesis [12, 28] based 
on Gaussian spatial covariance functions, is remarkably better 
than the SI reconstruction despite the facts that the LES fields 
being reconstructed did not actually have a Gaussian spatial 
covariance and may not exactly satisfy the assumptions of this 
hypothesis. 
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Appendix A. Elements of matrices R md and R dd 


where i = 1,2,...,/,/ = 1,2,..., 3 J, k = 1, 2, ..., Q, 
r e L t , B tt , B uu , B vv and B uv are the spatial-temporal 
covariance functions of the corresponding fields labeled by 
the subscripts, and the r 7 are the spatial points within the 
tomographic plane at which the fields are being reconstructed; 
these points remain fixed during the integration. 

Similarly, an expression for the covariance matrix (dodj) 
between the noise-free data at time t n and time t k can be 
calculated: 


-- f d if d/'j 4^ 

JLi Jl p 147b 


Co (4)Co (t n ) 
0 tk)To(t n ) 


(4)*70/? (4)) — 

+ B uu (y, t n , r , t^SfxSp X + B vv (y, t n ,Y , tfc)SiySpy 
+ B uv (y, t n ,Y, t^SfxSpy + B vu (y, t n ,Y , tji)SiySp X 


B TT {r, 4;r',4) 


(A.4) 


where i, p = 1, 2,...,/, n, k = 1, 2, ..., Q, r e L t , and 
y' e L p . Note that B vu (y, 4; y' , t 2 ) = B uv (y' , t 2 ; r, 4), and 
similarly for other fields. When deriving equations (A.3) and 
(A.4), it is assumed that B Tu = B Tv = 0. 


In this appendix, formulae for R dd and R md matrices are 
presented. 

The data vector d, given by equation (20), contains noise 
due to errors in the travel-time measurements, transducer 
positions and the reconstructed mean fields. Therefore, the 
data vector d can be presented as d = do + Ad, where 
d 0 = [qo(4); qo(4); • • •; qo(^e)] is a noise-free data vector, 
whose elements are given by equation (18), and Ad = 
[|(4); |(4); ...; £(62)] * s a vector of uncertainties with zero 
mathematical expectation: (Ad) = 0. Assuming that Ad is 
independent of m and do, one has 

Rmd = <md r > = (mdp) + (oiKAd 7 } = (mdj), (A.l) 

Rdd = <dd r ) = (dodj) + (AdAd r ). (A.2) 

As one can see, noise in the data does not affect the model-data 
noise-free covariance matrix (R mdo = R m d) but changes the 
data-data noise-free matrix R dodo . Namely, one should know 
the noise covariance matrix = (AdAd r ). The elements 
of this matrix are calculated in section 4.2. 

To calculate the elements of noise-free matrices (mdj) 
and (dodj) in equations (A.l) and (A.2), one should take into 
account equation (18): 

(mj(t 0 )q 0i (t k ))= [ dl\^L {mj(to)Tir j k)) 

JLi [ 2/o(4) 

+ (mj(t 0 )u(Y, t k ))s ix + (mj(t 0 )v(Y, t k ))s iy 

f Ll dl \_M^) BTT ^ Tj,to; r ’ tk ^\ ’ if 1 ^ J ^ -A 

Ili j-j^ to, r, t k )six + B uv (y j—j, to, r, t k )siy], 

= if J + 1 ^ j ^ 2J, 

I L{ j—2J ■> to-> U t k ^)Six + B vv (y j _ 27 , to, T, t k ^S(y], 

if 2J + 1 ^ j ^ 3 J, 


Appendix B. Spatial-temporal covariance functions 


In this appendix, analytical formulae for spatial-temporal 
covariance functions of temperature and wind-velocity 
fluctuations are presented. These formulae were derived with 
the use of the locally frozen turbulence hypothesis [12, 28] 
which is a generalization of the widely used rigidly frozen 
turbulence hypothesis (Taylor’s hypothesis). In the latter it is 
assumed that each point of a turbulent eddy is advected with a 
constant velocity. As a result, the eddy remains ‘frozen’ as it 
moves. In contrast, in the locally frozen turbulence hypothesis, 
it is assumed that each point of the eddy can move with its own 
velocity. Therefore, the turbulence is no longer ‘frozen’ since 
turbulent eddies can arbitrarily change their shape. If the wind 
velocity was constant, the locally frozen turbulence hypothesis 
coincides with Taylor’s hypothesis. 

For the two-dimensional case, the formulae for the 
spatial-temporal covariance functions of temperature and 
wind-velocity fluctuations at points r ' and y" and times t' and 
t" take the following form [12]: 


Btt(p ,t) = a? exp 


(p-Vo( f')T) 2 
l 2 T 


Buu(p,r) = d v exp 

(Py - V 0 (t')z) 2 


(p-Vo(f')T) 


'W\ 2 ”l 




X 1 - 


l 2 v 


(p - Vo(r')r) 




7 2 

l v 


B uv {p,z) = o v exp 

iPx ~ uo(t')r)(p y - Do (t')r) 

x 72 

L v 


(B.l) 


(B.2) 


(B.3) 


(A.3) 
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where p = r"—r', r = t" — t', Vo is a spatially averaged wind 
velocity with components (mq, vq), and the effective variances 
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dj and dy 

, and the squares of correlation lengths l\ 

and ly are 

[12] 

given by 






G 2 ~ 

2 

Gj 


Ij = l\ + 2 (TyT 2 , 

. (B.4) 

[13] 

D j 1 — 


3/2 ’ 

& v = 

”5 


72 i2 , 0 2_2 

V V + 2<7yT 

. (B.5) 

[14] 

[i+2(^) 2 : 

j 3/2 ’ 


Here, a T , 

Gy , Ij and ly 

are variances and correlation lengths 

[15] 

of the fluctuations at the same time (r =0). Formulae for B uv 

[16] 


and B vv are symmetric to those given by equations (B.2) and 
(B.3). Note that the space variations of the full wind velocity 
field V(r, t) in the locally frozen turbulence hypothesis are 
reflected in formulae (B.1)-(B.5) through Gy. 

In the limiting case g v t/ min{l v , l T j 0 and time- 
independent Vo, these formulae yield the results for rigidly 
frozen turbulence [11]. It follows from equations (B.l)- 
(B.5) that the dependence of the spatial-temporal covariance 
functions on r in the locally frozen turbulence hypothesis 
is manifested in three effects: the spatial arguments of the 
covariance functions are shifted by the vector the 

effective variances of the fluctuations and Gy decrease, 


and the effective correlation lengths l T 
increases. 
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